rm(list = ls())

source('param_general.R')
source('./SwitchCosts/stata_results_switch.R')


load(paste0(datDir,'estimDat_minPpl_all_',minZipHospPpl,'.RData', sep = '')) 

## Compute zip code level hospital quality: only use first births
driveTimes <- estimDat[,list(time_current = min(time_current),time_current_check = max(time_current)),list(hosp_id_alt, pat_zip)]
estimDatUse <- estimDat[age06 <= ageCut06 & birthNum ==1, ]

grpShrDat <- estimDatUse[,list(nGrpHosp = sum(obs_choice)), by = list(pat_zip, hosp_id_alt)]
grpShrDat[,nGrp := sum(nGrpHosp), by = list(pat_zip)]
grpShrDat[,hospShare := nGrpHosp/nGrp]


grpShrDat[,nGrp := sum(nGrpHosp), by = list(pat_zip)]
grpShrDatReshape <- merge(subset(grpShrDat, hosp_id_alt != 0),
                          subset(grpShrDat, hosp_id_alt == 0,
                                 c('pat_zip', 'hospShare')),
                          by = c('pat_zip'))

setnames(grpShrDatReshape, c('hospShare.x','hospShare.y'), c('hospShare','hospShareOut'))
grpShrDatReshape[,deltaHat := log(hospShare) - log(hospShareOut)]

nrow(driveTimes)
setkey(driveTimes,hosp_id_alt, pat_zip)
driveTimes <- unique(driveTimes)
nrow(driveTimes)
grpShrDatReshape <- merge(grpShrDatReshape, driveTimes, by = c('hosp_id_alt', 'pat_zip'))

grpShrDatReshape <- grpShrDatReshape[nGrp >= xiSmallGrp,]
grpShrDatXi <- data.table(ldply(coefDistList,
                                function(coef) copy(grpShrDatReshape)[,xiHat := deltaHat - coef * (time_current)]))
grpShrDatXi <- grpShrDatXi[is.finite(xiHat) == TRUE,]

##pull in aggregate hospital quality
aggDeltDatList <- list()
i <- 1
for (rt in c('lagdv','stlogit','fe')){
    for (ps in c('honoreSub','baseSub')){
        aggDeltDatRaw <- fread(paste0(datDir,'hosp_dumvar_agg_',rt,'_',ps,'.csv'))
        aggDeltDatList[[i]] <- aggDeltDatRaw[V2 %like% 'hosp_id_alt',list(hosp_id_alt = str_extract(V2,"[1-9][0-9]*"),aggXi = V3,popSet = ps,regTyp = rt)]
        i <- i+1
    }
}

aggDeltDat <- rbindlist(aggDeltDatList)[!is.na(hosp_id_alt)]
aggDeltDatWide <- dcast.data.table(aggDeltDat,hosp_id_alt~popSet+regTyp,value.var = 'aggXi')
cor(aggDeltDatWide[,c(2,3,4),with=FALSE])

save(grpShrDatXi, aggDeltDat,file = paste0(datDir,'deltaDat_minPpl_all_',minZipHospPpl,'.RData'))
